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Abstract. Public transport routes sharing the same grid of streets and tracks are often found to proceed 
in parallel along shorter or longer sequences of stations. Similar phenomena are observed in other 
networks built with space consuming links such as cables, vessels, pipes, neurons, etc. In the case of 
public transport networks (PTNs) this behavior may be easily worked out on the basis of sequences 
of stations serviced by each route. To quantify this behavior we use the recently introduced notion 
of network harness. It is described by the harness distribution P(r,s): the number of sequences of s 
consecutive stations that are serviced by r parallel routes. For certain PTNs that we have analyzed we 
observe that the harness distribution may be described by power laws. These power laws observed 
indicate a certain level of organization and planning which may be driven by the need to minimize the 
costs of infrastructure and secondly by the fact that points of interest tend to be clustered in certain 
locations of a city. This effect may be seen as a result of the strong interdependence of the evolutions 
of both the city and its PTN. 

To further investigate the significance of the empirical results we have studied one- and two- 
dimensional models of randomly placed routes modeled by different types of walks. While in one 
dimension an analytic treatment was successful, the two dimensional case was studied by simulations 
showing that the empirical results for real PTNs deviate significantly from those expected for randomly 
placed routes. 
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INTRODUCTION 



A variety of different phenomena have in recent years been analyzed in the context of complex 
network theory [1, 2]. Usually, the focus is on the network topology while the study of specific 
features as e.g. the network load or real-space correlations is mostly left asidej^] 
Analyzing statistical properties of public transport networks (PTNs) [3, 4, 5] a so-called harness 
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effect which results from such spatial correlations has recently been proposed [5, 6]. The latter 
may be observed for networks on which a set of walks or paths is defined. Indeed, public transport 
routes sharing the same grid of streets and tracks define such a set of walks on this grid. Often these 
are found to proceed in parallel along shorter or longer sequences of stations. Similar phenomena 
are observed in other networks built with space consuming links such as cables [7], vessels [8], 
pipes [9], neurons [10], etc. In the case of PTNs, a quantitative description of sequences of stations 
that are served by several routes may be performed in the form of the harness distribution P(r,s): 
the number of maximal found on the network sequences of s consecutive stations that are serviced 
by r parallel routes. In a recent empirical analysis of the harness effect on PTNs of 14 major cities 
in the world [5, 6] it was found that for certain cities the behavior of the function P(r, s) is described 
by power laws. These observed power laws indicate a certain level of organization and planning 
which may be driven by the need to minimize the costs of infrastructure and secondly by the fact 
that points of interest tend to be clustered in certain locations of a city. 

In the present paper, to further investigate the significance of the empirical results we have 
studied one- and two-dimensional models of randomly placed routes modeled by different types 
of walks. The setup of the paper is the following. First we recall some results of the empirical 
analysis. We then proceed solving a simple Id model of a growing non-correlated network that 
allows an analytical treatment. Finally we analyze different types of growing networks on the 
2d square lattice and perform numerical simulations and measure the resulting harness effect as 
function of the model parameters. Conclusions are given in the last section. 

RESULTS OF AN EMPIRICAL ANALYSIS 

In a previous study we have analyzed the PTNs of 14 major cities that have different geographical, 
historical and cultural background [5, 6]. For some of these we observed non- vanishing harness 
distributions P(r,s) > even for long sequences s and high numbers of routes r. This is what we 
call a "strong" harness effect (examples are Sao Paolo, Hong Kong, Istanbul, Los Angeles, Rome, 
Sydney, Taipei, Moscow, London; some are shown in Fig. 1). For other PTNs the maximal values 
of s and r with P(r,s) > were found to be smaller than 10 (Berlin, Paris, Dallas, Duesseldorf, 
Hamburg) - this we call a "weak" harness effect (Fig. 2). It is important to note that the division 
into these two classes does not correlate with either the average number of routes R in the PTN 
or their average length S. Another result is that similar to the node-degree distributions [5, 6] we 
observe that the harness distribution P(r, s) for some of the cities (Sao Paolo, Hong-Kong, Istanbul, 
Los Angeles, Rome, Sydney) may be described by a power law: 

P(r,s)~r~ Ys , for fixed s, (1) 

whereas the PTNs of other cities (Taipei, Moscow, London, Hamburg) are better described by an 
exponential decay: 

P{r,s)~e~ r/fs , for fixed s. (2) 

We illustrate this behaviour in Figs. la,b showing the harness distribution for Istanbul and for 
Taipei. In some cases (e.g. for Rome and Los Angeles) there is a crossover between regimes (1) and 
(2) at larger s as shown this for the PTN of Los Angeles (Fig.lc). Here, one can see that for small 
values of s the results are better described by a power law (1). With increasing s a tendency to an 
exponential decay (2) appears (Fig. Id). This is less obvious for other cities analyzed, however in 
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FIGURE 1. s-cumulative harness distribution P c (r, s) as function of r for fixed s. Log-log for Istanbul (a) 
and Taipei (b) s = 2,5, 8, 1 1, 14. Log-log for Los Angeles s = 2,4,6,8, 10 (c) and log-lin for Los Angeles 
£ = 22,24,26,28,30 (d). 



all cases the harness distribution P(r,s) as function of r decays faster for longer sequence lengths 
and while also attaining a more pronounced curvature. 

Note that in Fig. 1 we plot ^-cumulative distributions P c (r,s) where a sequence with maximal 
length s = 9 will be counted once as a sequence of length s = 9 and twice as a sequence of length 
5 = 8 etc: 

5 

P c (r,s) = Y,(s+l-s)P(r,s) (3) 

s=s 

It may be surprising that these curves e.g. for Taipei (Fig. lb) intersect for low values of r. We will 
discuss this effect below. 

For PTNs for which the harness distribution follows a power law (1) the corresponding exponents 
y s are found in the range of % = 2 4. For those distributions with an exponential decay the scale 
f s (see eq.(2) varies in a range f s = 1.5 4-4. The power laws observed for the behavior of P(r,s) 
indicate a certain level of organization and planning which may be driven by the need to minimize 
the costs of infrastructure and secondly by the fact that points of interest tend to be clustered in 
certain locations of a city. Note that this effect may be seen as a result of the strong interdependence 
of the evolutions of both the city and its PTN. We want to emphasize that the harness effect is a 
feature of the network given in terms of its routes but it is invisible in any of the complex network 
representations of public transport networks presented so far, such as L-space [3], P-space [11] or 



B-space [5, 12]. It is possible, that the notion of harness may be useful also for the description of 
other networks with similar properties. On the one hand, the harness distribution is closely related 
to distributions of flow and load on the network. On the other hand, in the situation of space- 
consuming links (such as tracks, cables, neurons, pipes, vessels) the information about the harness 
behavior may be important with respect to the spatial optimization of networks. A generalization 
may be readily formulated to account for real- world networks in which links (such as cables) are 
organized in parallel over a certain spatial distance. While for the PTN this distance is simply 
measured by the length of a sequence of stations, a more general measure would be the length of 
the contour along which these links proceed in parallel. 

For the cities observed no correlation appears to occur between the harness distribution behavior 
and other well-known network characteristics that were analyzed, as for example the node-degree 
distribution of PTNs [5]. 
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FIGURE 2. ^-cumulative harness distribution P c (r,s) as function of r for fixed s for Paris (s = 2 + 8). 
P c (r, s ) = if s > 2 for r > 7 and if s > 7 even for r > 2. 

However, the extent to which harness properties are expressed may obviously play a role for the 
attack vulnerability of a PTN. Interestingly, our previous investigations have shown that the Paris 
PTN is most resilient to any type of random or directed attacks (in terms of percolation concepts) 
among all analyzed PTNs [13]. At the same time it exhibits the "weakest" behavior with respect 
to the harness effect (Fig. 2). One may expect such a result: routes that do not share the same 
streets are more resilient. However, for other cities no apparent correlation between harness effect 
behavior and their vulnerability has been found so far. 

As noted above the interesting question to answer is: are there any structural evolution purposes 
behind this effect, or can it be found just as well within simple random scenarios. 

ANALYTIC RESULTS AND MODELING IN ID 

Let us first investigate a network model with routes placed randomly in one dimensional space. 
Although being very simple this model can mimic a harness effect and as we will see below, it 
allows an analytical solution. The model is formulated in the following way: 

The left terminals of R routes of length S are placed at random on a line with N sites, with periodic 
boundary conditions. E.g. in Fig. 3 we show two routes of length 5=15 with left terminals at x = 
and x = — 8. We define the route density as p = R/N routes per site. 
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FIGURE 3. R = 2 routes given as simple sequences of S = 15 consecutive sites are placed at random on 
a line with N sites. 



The distribution of left-terminals on a given site, e.g. site x will be P x (r): the probability that r 
routes have their left-terminal on site x 

P ' W= (r)4 )r ' (1 -^ )B " r ' (4) 

where the first term counts the number of ways to select r of R routes, the second term is the 
probability that the r left terminals lie on site x and the third term is the probability that no left 
terminal of an unselected route lies on site x. In other words, by definition P x (r) is a binomial 
distribution. For N — > °°, but fixed p —R/N this distribution has the limiting behavior of a Poisson 
distribution: 

P x (r)^e~P-^. (5) 

Let us now calculate the probability that there is a sequence of maximal length s and maximal 
width r of r routes in parallel between sites x = and x = s. This implies that at least one of the r 
routes starts at x = and at least one of the routes ends at x = s — 1. The latter route then starts at 
x = — (S — s) = —5 (Fig. 3). The other r — 2 routes may start anywhere in between — s < x < 0. We 
denote the number of routes starting at x = as ro, those starting at x < as r_ x . 

In the limit R^> r and W> jwe may consider Po(ro), Pi(r\), ... ,Ps(ri) as independent proba- 
bilities (this is not true for small systems, however, as we see later the correlation between these 
probabilities is negligible for the cases studied here. The overall probability to find a sequence of 
length s and width r starting at x = is then the sum over all combinations leading to the result: 
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PoM = £ P(r ).P(r 1 )....-P(r,)= £ e~' s P 
{r f },fiD>l,rr>l {n},ro>\,r s >\ 
(r) i 
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where L[^| denotes a sum over {r ; } with r = tq + r\ + ... + r^. 
Now, without the conditions ro > 1 and r$ > 1 this sum can be derived from: 

(s - +ir = (1 + 1 + 1+ ... +ir = |_^L_. (7) 

The sum with these conditions however can be written as: 

to i to j to x 
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FIGURE 4. s-cumulative harness distribution P c (r,s) as function of r for fixed s. Log-lin scale, a) 
Comparing the analytical solution and numerical simulations for N = 10000, R = 1000, S=10; b) for 
analytical solution for N = 10000, R = 2000, S=20; c,d) comparing the analytical solution with empirical 
results for the Moscow PTN (R an = Rmoscow = 679, S an = Smoscow = 22, N = 5250) for different s normalized 
by P(2,s). 



Thus: 

P ( r , s ) = g (-*P) .p'. [(s+ l) r -2(s) r +(s- \) r ]/r\. (9) 

With this formula we count sequences that start at x = 0. To receive the overall probability this 
is to be multiplied by N. 

P(r,s) = N ■ e { -^ ■ p r ■ [(s + l) r - 2{s) r + (s- \) r ]/r\. (10) 

Simple arithmetic allows us to calculate the ^-cumulative distributions P c (r,s), see eq. (3). 

As mentioned above, the probabilities we use are appropriate for infinite systems. To test their 
validity for finite cases we performed some simple simulations. It turns out that the average results 
fit this formula with very good accuracy even for small N, for example N = 10 and of course for 
any larger values of N (Fig. 4a). For a large range of parameters R, S, N the behavior of P(r,s) 
looks similar to what is shown in Fig. 4b. For high overall density p • S we observe that curves 
for different s intersect at small r as for some real- world PTNs. In one dimension, this effect has 
the following explanation: when the space is overcrowded one will in general find more than two 
routes to overlap for small sequences of stations. 



This model has three parameters: the number of routes R, the route length S and the number of 
sites N. As is obvious from eq. (10) the harness distribution P(r,s) for all r and s is close to a 
Poisson decay (5) for any set of parameters. Therefore PTNs with an exponential behavior of the 
harness distribution P(r, s) may be compared with the results of the one-dimensional approach. As 
an example we compare the normalized harness distribution for Moscow and the one-dimensional 
set of lines (Fig. 4c,d), where the number of routes R and the route length S were chosen to match 
those of the Moscow PTN (5 set to the average route length). 

However, the quantitative results for all observed PTNs are several orders of magnitude higher 
than the result obtained with (9) in one dimension for the same R and S (for any N). Furthermore, 
PTNs that show a power law behavior (1) are even qualitatively different from the random ID 
approach. Another difference is that the harness distribution curves for different s are very similar 
in shape and both slope and curvature vary much less than for the PTN harness distributions. 

In the following we test a two-dimensional model with the simple simulations. 



SIMPLE 2D MODELING 



It is obvious that simply throwing at random lines parallel to the axis' of a 2d square lattice with 
periodic boundary conditions will lead to the original Id problem: If the lattice has X x X sites 
one would get 2X independent one-dimensional systems. However, it is not a priory clear what 
results one will find for more general sets of walks on a 2d square lattice. To work this out, 
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FIGURE 5. s-cumulative harness distribution P c (r, s) as function of r for fixed s (s = 2 8), for R = 500, 
5 = 30, X = 50. Routes are generated as: a)RW b)NRRW c)S AW. P(r, s) = if s > 2 for r > 1 1 and if s > 7 
even for r > 2. 



we implemented the following simulations. We work on a 2d X x X square lattice with periodic 
boundary conditions. On this lattice we chose a set of R walks each of length S (number of steps 
plus 1). The routes are built either as random walks (RW), non-reversal random walks (NRRW) 
that cannot reverse the previous step, or self-avoiding random walks (SAW), that may not intersect 
themselves. 

These models have three parameters: the number of routes R, the route length S and the lattice 
size X. We choose the first two parameters to match those of different real PTNs. 

Postponing a more detailed analysis to a separate publication we here summarize some of the 
main features of the harness distribution P c (r,s) of these models. Besides the finding that the 
harness effect is "weak", some similarity between the the harness effects seen in the three models 
is observed (Fig 5a,b,c). Curvature and slope evolve in a similar way. Also intersections between 
the curves for different s are found to occur at lower values of r in all cases. Differences are that the 
RW-generated networks demonstrate a "weaker" harness effect, while NRRW- and SAW-generated 
networks result in harness distributions P c (r,s) of similar order of magnitude. 

It turns out that for fixed R and 5, increasing the lattice size X, the harness distributions P(r,s), 
for all fixed r < R and s < S show non-monotonous behavior (Fig. 6). As function of X it first 
increases, and then after reaching a maximum it starts to decrease. Comparing with the empirical 
values found for real PTNs we observe, that for some PTNs the empirical values even for small r 
and s are significantly larger than the maximum that could be obtained with RW. For NRRW or 
SAW the empirical values are within the observed range, however only for a small interval of X. 
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FIGURE 6. Situations observed for the r- and ^-cumulative harness distribution P cc (f,s) at r = 2, s = 2 
as function of lattice size X. For all three of models with R = RHamburg, S = S Hamburg- The empirical value 
of P cc (2,2) for Hamburg is shown for comparison. 

This proves numerically the not surprising observation, that with any of the proposed random 
or quasi-random walks only a very "weak" harness effect may be obtained. In turn, this strongly 
indicates that for most of the observed cities the harness effect must have a structural background, 
that is not to be modeled by any of the random approaches taken here. 

Let us therefore analyze the harness effect for a model that intrinsically takes correlations into 
account using some ideas of non-equilibrium network evolution. [5]. 



AN EVOLUTIONARY MODEL IN 2D 



From previous studies of PTNs [4, 5] the following facts are relevant for our further discussion of 
PTN network evolution. The final network often displays features (a power law degree distribution) 




FIGURE 7. Simulated PTN maps. X = 300, R = 1024, 5 = 64. Parameter a = 0. a) b = 0. 1 b) b = 0.3 a) 
ft = 0.5. One can see how parameter b influences the spread of the network. 



that are compatible with a preferential attachment scenario [1]. Furthermore, the fractal dimension 
of public transport routes appears to match that of 2d self-avoiding walks [5], This leads to the 
following model. 

Given a 2d square lattice with side X and periodic boundary conditions, R routes each of length 
S will be created in the following way. A first route is built as self-avoiding walk. Let k x be the 
number of visits to lattice site x. The R —\ subsequent routes are constructed as SAWs with the 
following preferential attachment rules: 

a) choose a terminal station at xq with probability 

p~k XQ +a/X 2 ; (11) 




FIGURE 8. Simulated PTN maps. X = 300, R = 1024, S = 64. Parameter b = 0.5. a) a = 15 b) a = 20 a) 
a = 500. One can see how parameter a influences the number of clusters. 

b) choose any subsequent station x of the route with probability 

p~k x + b. (12) 

c) If route intersects itself, discard it and return to step a). 
Repeat steps a)-c) until R routes are created. 

We call this model the model of mutually interacting self-avoiding walks (MI SAW). It imple- 
ments a preferential attachment scenario for networks that are built from chains (routes in our 
case). On the figures presented (Fig. 7-8) one can see how the parameters a and b influence the 
distribution of the routes on the lattice. Parameter a controls the number of disconnected clusters, 
while parameter b is related to the spread of each cluster. If both parameters equal zero, all R routes 
are restricted to the sites occupied by the first route. 

To present the obtained results we show the harness distribution P c (r,s) for different values for 
the model parameters a and b (Fig. 9). It turns out that for choosing different values we may 



1000 



-400 



10 




(a) 



1 000 



100 



10 




1 000 



100 



10 




(b) 



1 000 



100 



T5 20^ ^5 ' 




5 Wb; #" h ' ^) DD 25 



30 



(c) (d) 

FIGURE 9. s-cumulative harness distribution P c (r, s) as function of r for fixed s (s = 2,4,6,8, 10), for MI 
SAW R = 500, S = 30, X = 50. Log-log for a) a=0, b=0.5 and b) a=500, b=0.1. Log-lin for c) a=15, b=0.5 
and d) a=500, b=0.5. 



either find a behavior that corresponds to a power law distribution (Fig. 9a,b) or a behavior that 
corresponds to an exponential decay (Fig. 9c,d). 

Also other features of real- world PTNs may be reproduced. With increasing the sequence length 
s the harness distribution P(r, s) as a function of the number of routes r attains a stronger curvature 
and a steeper slope. Thus, with increasing sequence length s there is a cross-over from the power 
law regime (1) to an exponential regime (2). 

In the same way as observed for real-world PTNs in the Id approach and also for simple 2d 
models on the square lattice the ^-cumulative harness distributions P c (r, s) as function of the number 
of routes r intersect at low values of r. In one dimension we observed this behaviour because of an 
overcrowding effect, see above. This also occurs to some extent in the 2d case, however, on the 2d 
lattice there is an additional combinatorial effect that leads to this result: the number of different 
possible configurations of sequences with length 2 is smaller than the corresponding number of 
different possible configurations with length 3. To summarize one can say that the MI SAW model 
reproduces a large amount of the empirically observed behavior of harness distributions of PTNs. 



CONCLUSIONS 



Harness phenomena may be observed in different networks built with space consuming links such 
as cables, vessels, pipes, neurons, etc. The present analysis may possibly be applied also to such 
types of networks. In the particular case of the PTNs that we have analyzed we observe that in 
some cases the harness distribution may be described by power laws. These observed power laws 
indicate a certain level of organization and planning which may be driven by the need to minimize 
the costs of infrastructure and secondly by the fact that points of interest tend to be clustered in 
certain locations of a city. This effect may be seen as a result of the strong interdependence of the 
evolutions of both the city and its PTN. 

To further investigate the significance of the empirical results we have studied one- and two- 
dimensional models of randomly placed routes modeled by different types of walks. While in one 
dimension an analytic treatment was successful, the two dimensional case was studied by extensive 
simulations. 

Our main results are the following: 

• A one dimensional model for harness distributions was solved analytically. 

• Exponentially decaying harness distributions may be reproduced by the Id approach. 

• Simple random placement of RW, SAW or NRRW on a two dimensional square lattice result 
in weak harness distributions; in the RW case much weaker than for real PTNs. 

• The s-cummulative distributions for different s intersect at low values of r for all models due 
to combinatorial reasons. 

• A model of mutually interacting SAWs reproduces many of the empirically observed features 
of harness distributions. 
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